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ABSTRACT 


A matheinati cal  model  of  the  local  transient  diffusion-controlled  reaction 
between  initially  unmixed  species  is  presented.  It  is  intended  ultimately  as 
a computational  "molecule"  to  be  imbedded  in  direct  simulations  of  larger 
scale  reacting  flows.  The  model  consists  of  an  interacting  three-dimensional 
strain  vortex  field  which  exactly  satisfies  the  Navi er-Stokes  equation,  an 
analytically  determined  Lagrangian  representati on  of  the  mixing  process  and 
convection-diffusion  equations  for  the  reacting  species  in 
Lagrangian  coordinates.  The  length  scale  established  by  the  stretching  of  the 
vorticity  field  is  shown  to  be  directly  relatable  to  the  Kolmogoroff  scale  if 
the  local  strain  rate  has  a scale  consistent  with  laminar  boundary  layer 
theory  coordinates.  Results  are  shown  for  the  flow  pattern  and  the  induced 
mixing.  An  analytical  solution  to  the  convection-diffusion  equation  governing 
the  diffusion-controlled  reaction  is  derived.  The  solution  is  valid  for  large 
Schmidt  number  and  describes  the  evolution  of  any  initially  two-dimensional 
configuration  of  reactants.  A special  two-dimensional  case  of  this  model,  in 
which  vortex  strain  is  excluded  and  fuel  and  oxidizer  initially  occupy 
adjacent  half-spaces,  is  also  analyzed.  This  problem  was  originally 
formulated  by  F.  Marble,  who  treated  it  in  a very  different  manner  from  that 
described  here. 


INTRODUCTION 


The  numerical  simulation  of  turbulent  reacting  flows  has  become  increasingly 
viable  in  recent  years.  More  powerful  computers,  advances  in  computational 
techniques,  and  clearer  physical  insight,  have  combined  to  make  possible  the 
calculation  of  at  least  some  aspects  of  turbulent  combustion  directly  from  the 
governing  equations. 

Despite  this  progress,  many  serious  problems  remain  unsolved.  Although  large 
scale  features  of  turbulent  flows  seem  amenable  to  computation,  the  prediction 
of  local  reactant  consumption  and  energy  release  rates  requires  a description 
of  the  flow  at  length  and  time  scales  sufficiently  small  for  molecular 
diffusion  to  be  effective.  These  scales  are  well  below  the  resolution  limit 
of  any  finite  difference  simulation;  particularly  since  the  inherently  three- 
dimensional  mechanism  of  vortex  stretching  must  be  accounted  for  in  the 
description  of  phenomena  at  this  level  of  detail.  In  fact,  it  is  the  vortex 
stretching  mechanism  that  is  responsible  for  the  fact  that  the  ratio  of 
dissipation  scale  to  macroscopic  length  is  proportional  to  (Re)~3/4  -j p fully 
developed  turbulent  flow.  Here  Re  is  a Reynolds  number  based  on  macroscopic 
velocity  and  length  scales.  The  description  of  flow,  mixing,  and  diffusion  of 
reactants  at  this  scale  is  the  subject  of  the  present  research. 

The  mathematical  model  presented  here  is  intended  ultimately  as  a 
computational  "molecule"  to  be  imbedded  in  direct  simulations  of  larger  scale 
reacting  flows.  In  this  spirit,  it  is  assumed  that  a three-di mensi onal 
time-dependent  calculation  capable  of  resolving  the  largest  length  and  time 


scales  of  the  relevant  phenomena,  is  being  performed.  The  resolution  limits 
of  this  calculation  are  too  coarse  to  directly  calculate  the  small  scale 
mixing  and  diffusion  which  actually  control  the  reaction.  The  large-scale 
simulation  does,  however,  define  the  environment  in  which  the  phenomena  of 
interest  evolve. 

Progress  made  this  year  has  been  reported  at  the  1985  AFOSR/ONR  Contractors' 
Meeting  on  Turbulent  Combustion  on  23-25  July  in  Pasadena,  California  in  a 
talk  entitled  "Time-Dependent  Simulation  of  Turbulent  Combustion"  by 
H.  R.  Baum  and  R.  G.  Rehm;  and  at  the  Eastern  States  Section/Combustion 
Institute  Fall  1985  Meeting  on  November  4-6,  1985  in  Philadelphia, 

Pennsylvania  in  a talk  entitled,  "Time-Dependent  Simulation  of  Small-Scale 
Turbulent  Mixing  and  Reaction",  by  H.  R.  Baum,  D.  M.  Corley  and  R.  G.  Rehm.  A 
paper  based  on  this  work  entitled,  "Time-Dependent  Simulation  of  Small-Scale 
Turbulent  Mixing  and  Reaction",  by  H.  R.  Baum,  D.  M.  Corley  and  R.  G.  Rehm, 
has  been  submitted  to  The  Twenty-First  Internati onal  Symposium  on  Combustion, 
which  is  to  be  held  in  Munich,  West  Germany  on  August  3-8,  1986.  A copy  of 
this  paper  appears  as  Section  II.  In  Section  III  additional 
analyti cal /numeri cal  research  on  the  special  two-dimensional  case  described  as 
"the  Marble  problem"  is  presented;  this  is  work  in  progress  and  will  be 
published  later  when  it  is  completed. 
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ABSTRACT 


A mathematical  model  of  the  local  transient  diffusion-controlled  reaction 
between  initially  unmixed  species  is  presented.  It  is  intended  ultimately  as 
a computational  "molecule"  to  be  imbedded  in  direct  simulations  of  larger 
scale  reacting  flows.  The  model  consists  of  an  interacting  three-dimensional 
strain  vortex  field  which  exactly  satisfies  the  Navier-Stokes  equation,  an 
analytically  determined  Lagrangian  representation  of  the  mixing  process  and 
convection-diffusion  equations  for  the  reacting  species  in 

Lagrangian  coordinates.  The  length  scale  established  by  the  stretching  of  the 
vorticity  field  is  shown  to  be  directly  relatable  to  the  Kolmogoroff  scale  if 
the  local  strain  rate  has  a scale  consistent  with  laminar  boundary  layer 
theory  coordinates.  Results  are  shown  for  the  flow  pattern  and  the  induced 
mixing.  An  analytical  solution  to  the  convection-diffusion  equation  governing 
the  diffusion-controlled  reaction  is  derived.  The  solution  is  valid  for  large 
Schmidt  number  and  describes  the  evolution  of  any  initially  two-dimensional 
configuration  of  reactants. 
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INTRODUCTION 


The  numerical  simulation  of  turbulent  reacting  flows  has  become  Increasingly 

viable  in  recent  years.  More  powerful  computers,  advances  in  computational 

techniques,  and  clearer  physical  insight,  have  combined  to  make  possible  the 

calculation  of  at  least  some  aspects  of  turbulent  combustion  directly  from  the 

governing  equations.  For  example,  finite  difference  simulations  of  inviscid 

turbulent  mixing  layers  in  two  dimensions  and  axially  symmetric  jets  have 

been  performed  at  fairly  high  resolution.  Less  detailed  calculations  of 

simple  transport  by  buoyant  plumes  in  enclosures  have  been  carried  out  using 

finite-difference  and  particle-tracking  methods  in  both  two  and  three 
3 4 

dimensions  » . Finally,  the  random  vortex  method  has  been  combined  with  a 
flame  front  algorithm  and  finite  difference  techniques  to  produce  elegant 
simulations  of  preialxed  turbulent  combustion  in  two  dimensions'*^. 

Despite  this  progress,  many  serious  problems  remain  unsolved.  Although  large 

scale  features  of  turbulent  flows  seem  amenable  to  computation,  the  prediction 

of  local  reactant  consumption  and  energy  release  rates  requires  a description 

of  the  flow  at  length  and  time  scales  sufficiently  small  for  molecular 

diffusion  to  be  effective.  These  scales  are  well  below  the  resolution  limit 

of  any  finite  difference  simulation;  particularly  since  the  inherently  three- 

dimensional  mechanism  of  vortex  stretching  must  be  accounted  for  in  the  / 

description  of  phenomena  at  this  level  of  detail.  In  fact,  it  is  the  vortex 

stretching  mechanism  that  is  responsible  for  the  fact  that  the  ratio  of 

dissipation  scale  to  macroscopic  length  is  proportional  to  (Re)~3/4  in  fully 

1 3 

developed  turbulent  flow  . Here  Re  is  a Reynolds  number  based  on  macroscopic 
velocity  and  length  scales.  The  description  of  flow,  mixing,  and  diffusion  of 
reactants  at  this  scale  is  the  subject  of  the  present  paper. 
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The  mathematical  model  presented  below  is  intended  ultimately  as  a 
computational  "molecule"  to  be  imbedded  in  direct  simulations  of  larger  scale 
reacting  flows.  In  this  spirit,  it  is  assumed  that  a three-dimensional 
time-dependent  calculation  capable  of  resolving  the  largest  length  and  time 
scales  of  the  relevant  phenomena,  is  being  performed.  The  resolution  limits 
of  this  calculation  are  too  coarse  to  directly  calculate  the  small  scale 
mixing  and  diffusion  which  actually  control  the  reaction.  The  large-scale 
simulation  does,  however,  define  the  environment  in  which  the  phenomena  of 
interest  evolve.  This  environment  can  be  defined  as  follows:  Classical 

analyses  of  the  kinematics  of  an  arbitrary  flow^  show  that  the  velocity  field 
in  the  neighborhood  of  any  point  can  be  decomposed  into  a uniform 
translational  velocity,  a rigid  body  rotation,  and  a pure  strain.  This 
decomposition  can  be  readily  performed  for  most  numerical  simulations. 
Similarly,  the  thermodynamic  state  expressed  in  terms  of  the  pressure  and 
conserved  scalar  quantities,  (e.g.  Schvab-Zeldovich  variables)  are  assumed 
known.  At  high  Reynolds  numbers,  the  conserved  scalars  are  primarily 
convected  at  computationally  resolveable  scales.  It  is  easy  to  track  the 
coordinates  of  a material  point  moving  with  the  computationally  resolveable 
fluid  motion  > . We  now  choose  a local  coordinate  system  in  this  frame  of 
reference. 

Viewed  from  the  moving  coordinate  system,  the  ambient  velocity  field  consists 
of  a linear  strain  field  and  a vorticity,  with  a specified  pressure.  There 
are  two  possibilities  for  the  conserved  scalar  fields.  If  the  nominal 
interface  between  reactants  is  macroscopically  determined  then  at  the 

resolution  limit  the  interface  is  taken  to  be  initially  planar.  If  It  Is  of 
subgrid  dimensions,  then  its  initial  shape  is  assumed  known.  Examples  of  both 
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types  are  shown  below.  The  fundamental  problem  can  now  be  defined  as  follows: 
Given  initial  conditions  for  the  "local  ambient  state"  (i.e.  strain  rate, 
vorticity,  pressure,  and  conserved  scalar  configuration),  determine  their 
subsequent  local  evolution  in  space  and  time. 

This  problem  is  broken  down  below  into  three  subsidiary  problems.  First,  the 
construction  of  the  local  velocity  field  produced  by  the  effect  of  the  strain 
on  the  vorticity  field.  Next,  the  calculation  of  a local  Lagrangian 
coordinate  system  so  that  small  scale  mixing  of  reactants  can  be  separated 
from  molecular  diffusion.  Finally,  solution  of  the  evolution  equation  for 
conserved  scalars  in  the  Lagrangian  coordinates.  In  each  case  simplifications 
in  geometry  and  thermophysical  properties  are  introduced  in  the  Interest  of 
maintaining  an  essentially  analytical  approach  yielding  readily  interpretable 
results.  The  consequences  of  these  simplifications  are  assessed  in  a 
concluding  section. 

A related  set  of  investigations  has  been  performed  by  Karagozlan  and 
3 9 

Marble  » . However,  both  the  general  approach  and  the  results  obtained  by 
them  are  significantly  different  from  those  reported  below.  More  detailed 
comments  are  included  in  subsequent  sections. 


THE  VELOCITY  FIELD 

Consider  a flow  in  which  the  vorticity  is  locally  aligned  in  the  x-direction. 
Fig.  1,  and  subjected  to  a uniform  strain  of  strength  a.  Then,  if  v^,  Vj-,  and 
V0  are  respectively  the  axial,  radial,  and  swirl  components  of  velocity,  the 
Navier  Stokes  equations  for  a constant-property  incompressible  fluid  have 
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solutions  of  the  form: 


Vx  = ax  ; Vf  = -ar/2 

9v0 
3t 

Equation  (1)  is  the  swirl  component  of  the  momentum  equation  for  the  fluid 
whose  kinematic  viscosity  is  v.  The  axial  vorticity  field  o)(r,t)  is  defined 
in  terms  of  V0(r,t)  by  the  expression 

m ^ (rv0) . (2) 

r 3r 

The  linearity  of  Eqs.  (1)  and  (2)  suggests  seeking  a solution  satisfying  the 
initial  condition 

u)(r,0)  = ro  J^6(r-ro)  (3) 

Here  6 denotes  the  Dirac  delta  function.  The  resulting  solution  for  V0  is  in 
effect  a Greens  function  for  Eq.  (1),  in  that  its  convolution  with  any  iniLiai 
axial  vorticity  distribution  generates  the  corresponding  solution  to  the 
Navier-Stokes  equations.  As  such,  it  contains  all  the  phenomena  implied  by 
the  equation;  i.e.  the  convection,  stretching,  and  diffusion  of  vorticity. 

To  proceed,  it  is  convenient  to  introduce  the  following  dimensionless 
variables : 

AAA 

V0  = ro  w(r,t) 

A 

r = ro  r (4) 

t = (a)“l  t 
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Then,  define  a new  dependent  variable  v,  and  independent  variables  p,T  in  the 
form 

A A 

p = r exp  (t/2) 

A 

T = e[exp(t)-l]  (5) 

A A 

w = exp(t/2)  v(p,t) 


In  terms  of  these  variables,  Eqs.  (l)-(3)  become: 

_n  = 1 ^ ( pax] 

9t  p 9p  3p 

X(p,0)  = 6(p-l)  (6) 

where 

X = 9 (pv) 

P 9p 

The  variable  p is  the  characteristic  coordinate  of  the  left  hand  side  of 
Eq.(l),  while  X is  a reduced  vorticity.  The  parameter  (e)“^  is  a local 
Reynolds  number  based  on  the  strain  rate  and  the  initial  radius  of  the  vortex 
tube  whose  evolution  is  described  by  X.  Note  that  the  first  of  Eqs. (6)  is  the 
heat  conduction  equation  with  cylindrical  symmetry.  The  solution  is  readily 
obtained  by  taking  Laplace  Transforms  in  The  result  is: 


X = (2t)  1 exp  [-(l+p2)/4T]  Io(p/2t) 


P 

V = (p)"l  / X (Po,t)  PodPo 
0 


(7) 


Here,  Iq  is  the  modified  Bessel  function  of  order  zero.  Equations  (5)  and  (7) 
constitute  the  desired  solution.  They  describe  the  evolution  of  an  initial 
shell  of  vorticity  concentrated  at  radius  ro  into  the  classical  steady  state 
stretched  vortex,  whose  velocity  is  given  by: 
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w(r,<»)  = (r)  1{  l-exp[-(r)2/4e  ] } 


(8) 


The  dimensionless  swirl  velocity  w defined  by  Eqs.  (4),  (5),  and  (7)  is 
plotted  in  Fig.  2 for  e = .002,  corresponding  to  a local  Reynolds  number  of 
500.  Note  that  the  stretching  mechanism  acts  to  steepen  the  velocity 
gradients  with  time  until  the  viscous  effects  become  strong  enough  to  provide 
the  steady-state  balance  described  by  Eq.  (8).  The  internal  length  scale  il  In 
these  solutions  is 


I = (v/a)l/2  (9) 

Thus,  (e)l/2  is  the  ratio  of  the  Internal  length  scale  to  rQ,  that  Imposed  by 
the  largest  scale  flow.  The  actual  magnitude  of  I depends  on  the  imposed 
strain  rate  ot.  It  is  interesting  to  observe  that  in  a laminar  boundary  or 
free  shear  layer  characterized  by  a macroscopic  length  L and  velocity  U,  the 
local  strain  rate  a can  be  characterized  by 

a ~ U(UL/v)1/2/L  (10) 

Inserting  the  estimate  given  by  Eq . (10)  in  Eq . (9),  the  Internal  scale  I can 
then  be  related  to  macroscopic  parameters  as: 

l/L  ~ (Re)"3/4 


(11) 


Re  = UL/v 

Equation  (11)  shows  that  under  circumstances  where  a fragment  of  laminar  shear 
layer  can  stretch  a vortex,  the  resulting  flow  will  quickly  collapse  to 
the  Kolmogoroff  or  dissipation  scale  Thus,  the  solution  described  here 

offers  a concrete  dynamic  realization  of  how  such  a small  scale  local  flow  can 
develop. 


The  fact  that  the  vorticlty  evolution  can  be  related  to  Lhe  cylindrical  heat 
conduction  equation  under  the  symmetry  conditions  implied  by  Eqs.  (1)  has  been 
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noted  by  several  authors  » > > . However,  when  solutions  investigated  are 

8 9 12 

restricted  to  Gaussian  forms  > » in  the  x,p  coordinate  system,  a rich 
variety  of  phenomena,  including  the  collapse  process  Itself,  is  lost.  This 
has  further  Implications  which  arise  in  the  study  of  mixing  and  diffusion, 
which  are  considered  next. 


REACTANT  MIXING 

Assume  that  the  fluid  consists  of  reacting  species  "A"  and  "B",  which  are 
initially  separated  on  either  the  macroscopic  scale  L or  the  local  scale  ro. 
For  simplicity  (and  consistency  with  the  flow  field)  the  reaction  is  taken  to 
be  Isothermal.  Let  n and  M be  the  stoichiometric  coefficient  and  molecular 
mass  of  each  species,  and  Y be  the  corresponding  mass  fraction.  Then,  if  Y(0) 
denotes  the  initial  unmixed  mass  fraction  of  each  species,  it  is  possible  to 
define  a normalized  Shvab-Zeldovich  variable  Z as  follows: 


Z = 2f(YA/nAMA)-(YB/nBMB)  - (1/2) [ (YA(0)/nAMA)  - YB(0)/nBMB]} 
•[(YA(0)/nAMA)  + 


(12) 


The  variable  Z satisfies  the  following  equation  and  initial  condition: 

az  + u*VZ  = DAZ  (13) 

at 

at  t = 0 Z = +1  (Species  A)  (14) 

Z = -1  (Species  B) 

Here  u is  the  velocity  vector  defined  in  Eq.  (1),  D is  the  diffusivity  of  each 
species,  and  A the  Laplacian  operator.  The  Initial  species  distributions  are 
assumed  to  be  independent  of  the  axial  coordinate  x.  This  restriction  is  not 
essential  to  the  analysis  which  follows,  and  it  greatly  simplifies  the  algebra 
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required  to  obtain  results  of  interest.  Then,  introducing  the  dimensionless 
variables  defined  by  Eq.  (4),  Eq,  (13)  can  be  written  in  the  form: 

3Z  r 3Z  w 3Z  el  32z  1 3Z  1 32z  I (15) 

A A A IaaAA  I 

3t  2 3r  a r 39  Sc  | 6r^  r 6r  r^  602  | 

Sc  = v/D 

The  variable  9 in  Eq.  (15)  is  the  azimuthal  coordinate,  increasing  in  the 

A 

direction  of  the  swirl  velocity  component  w (see  Fig.  1),  while  Sc  is  the 
Schmidt  number. 


The  left-hand  side  of  Eq.  (15)  can  be  thought  of  as  governing  the  small-scale 
mixing  of  the  reactants  by  the  local  velocity  field,  while  the  right-hand  side 
controls  the  true  molecular  diffusion.  It  is  desirable  to  cast  the  solution 
in  a form  which  illustrates  these  two  effects  explicitly.  This  can  be 
accomplished  by  introducing  Lagrangian  coordinates  x,p,(J).  The  variables  t and 
p are  defined  in  Eq.  (5),  while  (|)  is  defined  as: 

<j)  = 9 - (51/ae)  4'(p  ,t) 

T (16) 

iKp,t)  = (p)“l  / v(p,Tq)  dTo 
0 

Equation  (15)  in  Lagrangian  coordinates  becomes: 


3Z  , 3^Z  1 3Z  1 

— = (Sc)  1 { — - + - — + [-  + (Jl/ae) 
3t  2 p 3p  2 

3p  p 


i-r] 

3p 


- (Q/ae) 


8il;  3^Z  3 3ii;  3Z  1 3i|;  3Z 

[ — + — ( — — ) + - — — 1 } 

3p  3p3(j)  3p  3p  34>  P 3p  3(j> 


(17) 
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Note  that  in  this  system  of  coordinates  there  is  no  longer  any  convection,  as 
the  mixing  has  been  incorporated  in  the  transformation  from  Lagrangian  to 
Eulerian  coordinates.  Since  the  quantity  v(p,r)  is  given  analytically  by  Eqs. 
(7),  the  transformation  and  any  of  its  derivatives  can  be  evaluated  readily  to 
arbitrary  accuracy.  This  transformation  is  easily  generalized  to  three 
dimensions,  with  little  added  complexity  to  Eq.  (17). 

In  the  absence  of  diffusion  (Sc-»-«>),  Eq.  (17)  reduces  to  the  statement  that  the 
initial  configuration  of  the  reactants  does  not  change  in  Lagrangian 
coordinates.  Thus,  the  mixing  process  can  be  followed  by  mapping  each  point 
in  the  Eulerian  coordinates  at  a given  time  into  the  Lagrangian  space  and 
determining  which  species  occupied  that  point  originally.  This  has  been  done 
for  two  examples,  shown  in  Fig.  3 and  Fig.  4.  Figure  3 shows  the  mixing  of 
two  species  initially  occupying  separate  half  spaces  displaced  one  half  the 
initial  vortex  radius  from  the  vortex  center.  This  is  an  example  of  a 
macroscoplcally  resolveable  interface.  The  mixing  process  in  this  case  is 
seen  to  be  a relatively  large  scale  engulfment  followed  by  a continuous 
elongation  of  the  Interface  between  species.  Figure  4 shows  a "subgrid  scale" 
interface,  a circular  blob  of  one  reactant  displaced  one  half  the  vortex 
radius  with  respect  to  the  vortex  center.  This  time  the  mixing  seems  to  tear 
apart  the  blob,  with  two  highly  distorted  fragments  connected  by  a slender 
filament.  In  each  case,  the  local  Reynolds  number  (e)“^  was  taken  to  be  500, 
and  the  ratio  of  vorticity  to  strain  (.Q/a)  is  2.  While  the  Infinite  Schmidt 
number  limit  is  unrealistic  for  most  fluids,  it  is  an  excellent  approximation 
for  typical  smoke  aerosols. 
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MOLECULAR  DIFFUSION  EFFECTS 


Now  consider  the  solution  of  Eq.  (17)  for  finite  values  of  Sc.  In  the 
Lagrangian  representation,  the  vorticity  acts  to  enhance  the  diffusivity  in 
the  azimuthal  direction  in  Langranglan  coordinates.  The  degree  of  enhancement 
varies  greatly  with  p and  T.  Fig.  5 shows  the  variation  of  (J2/ae)i|;  with  r. 

A 

At  any  Instant  of  time,  however,  r and  p differ  only  by  a multiplicative 
constant,  so  each  individual  spatial  profile  has  the  same  shape  in  both 
systems  of  coordinates.  The  analogous  plot  in  Lagrangian  space  would  show 
each  profile  centered  at  p = l.  The  result  is  that  the  gradient  of  i|;  is  very 
large  in  an  annular  ring  centered  about  p = l.  Inside  this  ring,  ij;  vanishes 
exponentially,  while  outside  it  decays  algebraically  (Fig.  6). 


This  distribution  arises  from  the  fact  that  the  vorticity  is  confined  to  the 
annular  ring.  Thus,  there  is  no  angular  deformation  inside  the  annulus,  while 
the  deformation  outside  the  ring  is  irrotational . The  ring  thickens  with  time 
under  the  influence  of  viscosity,  until  the  inner  core  disappears.  Soon  after 
this  ocurs;  all  the  reactants  present  within  the  ring  (now  a disc)  will  be 
consumed . 


Since  the  coefficients  in  Eq.  (17)  are  independent  of  the  angular  variable  <> , 
it  is  useful  to  express  the  solution  in  terms  of  a Fourier  series. 

Z = I Znc(p  ,T)i;nc(P  cos  (n^) 


n 

+ I Zns(P  ,T)Cns(P  »T)sin  (n<J>). 
n 


(18) 


The  quantities  and  Z^s  sre  defined  to  be  the  Fourier  coefficients  that 

would  occur  in  the  absence  of  vorticity.  This  corresponds  to  the  solution  for 
diffusion  in  a uniform  strain  field,  a problem  which  has  already  received 
attentlon^'^ » ^5, 
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(19) 


/ dpipi  exp  [-Sc(p2+p2j^)/4x]  i^(Scppi/2t)| 


P 1 » 0) 
vZns  (P  1 > 0)/ 


Here,  are  the  Modified  Bessel  functions  of  order  n and  Znc(Pl>0)>  Zyis(pi,0) 
are  the  Fourier  coefficients  of  the  initial  distribution  of  reactants. 


To  proceed  further,  it  is  necessary  to  either  compute  the  remaining  factors 

?ns  numerically  or  make  an  additional  assumption  to  obtain  analytical 
results.  In  the  present  work,  the  latter  choice  is  made  and  it  is  assumed 
that  the  Schmidt  number  Sc  is  large.  Under  these  circumstances,  the  diffusion 
is  relatively  ineffective  in  changing  the  initial  distribution  except  in  the 

A 

annular  region.  Since  the  whole  vortex  collapse  process  occurs  in  times  t 
0(1),  then  T<<1  for  the  times  of  interest.  For  small  t,  it  can  be  shown  that 
^ has  the  form: 


>1^  = (p)"2  T 4'(n,T) 
n = (p-1)/2(t)1/2 


(20) 


The  quantity  'I  in  Eq.  (20)  is  a smooth  bounded  function  of  n which  tends  to 
unity  as  n->^  and  vanishes  exponentially  as 


For  convenience  in  what  follows  let  the  expansion  given  by  Eq.  (18)  be 

temporarily  recast  in  complex  exponential  form:  / 

Z = ^ exp  (in(}))Z„(p  ,T)?n(P  »t)  (21) 

n 
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Substitution  of  Eq.  (21)  into  Eq.  (17)  yields: 


9?n  9^?n  3^n 

Zn  = (Sc)-l  {Zn  + (2  

9t  3p2  6p 

3ip  9 Zn?n 

-in(S^/ae)[--  (—  {ZnCnl  + ) 

9p  9p  p 


Zn  9^n  9if; 

+ __  ) _ n2  (j^/ae)2  Z^Cn  ( " )^ 

P 9p  9p 

9 94^ 

+ --  (Zn?n  --  )]}  (22) 

9p  9p 


The  independent  variables  are  now  rescaled  as  follows: 

T = [ae(  Sc)  ^/2/q  ]t^* 

p-1  = [ae(Sc)l/2/s^]p* 
Equation  (20)  then  implies  the  scaling 
dJp 

— = q = [ae(Sc)^/2/j^]  1/2  q*(t*,p*) 

9p 


(23) 


(24) 


Since  the  reactants  are  initially  assumed  to  be  resolveable  on  a scale  of 
order  one  in  p,  then  Z^  and  its  derivatives  remain  0(1).  When  the  scaling 
laws  given  by  Eqs.  (23)  and  (24)  are  used  in  Eq.  (22),  the  result  Is: 

Zn  [ + n2(q*)2  = 0(Sc)“l/2  (25) 

9t* 


Thus,  the  leading  terms  in  the  expansion  for  large  Schmidt  number  can  be 
written  in  terms  of  p and  T as: 


T 94^ 

?nc  = ^ns  = exp  {-n2(Sc)“l  (Q/ae)^  f ( — )2  dx}  (26) 

0 9p 


Equations  (18),  (19),  and  (26)  then  constitute  the  solution  to  Eq . (17)  in  the 
annular  region.  They  are  in  fact  the  large  Schmidt  number  solution  every- 
where. To  see  this,  note  that  inside  the  annulus,  94^/9p  vanishes  exponen- 
tially and  Cnc  “ ^ns  ~ so  the  vortex  free  solution  is  recovered  exactly. 

For  values  of  p in  the  region  1 < p < <»,  Eq.  (22)  is  properly  scaled  In  p. 
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In  order  to  effect  a balance  for  times  t*  =0(1)  the  wave  number  n must  be 
rescaled  according  to 

n = [J2/ae(Sc)  ]n*  (27) 

This  is  merely  a Fourier  representation  of  a boundary  layer  in  the  angular 
variable  4>,  When  the  rescaling  is  performed,  the  same  terms  are  dominant. 

This  solution  can  be  interpreted  physically  as  saying  that  the  large  Schmidt 
number  prevents  significant  radial  diffusion  in  the  Lagrangian  coordinates. 
The  intense  vortlcity  sufficiently  enhances  the  azimuthal  diffusion  in  the 
annulus  so  that  only  relatively  low  wavenumber  angular  modes  are  required  to 
describe  the  solution.  In  the  irrotational  region,  the  diffusion  is  confined 
to  narrow  regions  near  the  interface  between  reactants  and  occurs 
preferentially  in  the  azimuthal  direction.  The  implications  for  the 
combustion  of  a fuel  blob  (see  Figure  (6))  are: 

(1)  The  fuel  in  the  annular  region  containing  vortlcity  is 
very  rapidly  consumed. 

(2)  By  the  time  the  vortex  collapses  to  its  final  steady  state,  the  only 
remaining  fuel  is  in  the  Irrotational  flow  region. 

(3)  Any  remaining  fuel  is  consumed  along  the  fuel-oxidizer  interface  in 
the  irrotational  region  over  a period  of  time  determined  by  the  extent  of 
the  fuel  blob  as  well  as  the  vortlcity  and  strain  rate. 
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CONCLUDING  REMARKS 


A mathematical  model  of  the  local  convection,  mixing,  and  diffusion  of 

initially  unmixed  reactants  in  a turbulent  flow  has  been  presented.  The  model 

is  based  directly  on  solutions  of  the  Navier-Stokes  equations  in  limiting 

geometrical  and  parametric  cases;  no  empirical  parameters  are  introduced. 

At  present,  the  model  cannot  be  regarded  as  complete.  In  particular,  more 

work  is  needed  to  extend  the  model  to  incorporate  three-dimensional  reactant 

configurations,  large  temperature  variations,  and  larger  molecular 

dif fusivities . It  is  currently  more  suited  to  the  analysis  of  acid-base 

reactions  in  turbulent  flow  than  to  combustion.  Still,  many  of  the  phenomena 

are  similar  in  both  reaction  systems,  and  the  study  of  acid-base  reactions  in 

jets  and  mixing  layers  is  an  ongoing  project  in  turbulent  combustion 
1 16 

research  > . Extended  versions  of  this  model,  employed  in  conjunction  with 

large-scale  flow  simulations,  offer  the  potential  for  predicting  turbulent 
reacting  flows  over  a much  wider  range  of  scales  than  is  possible  with 
supercomputers  using  conventional  numerical  methods.  This  statement  Is  likely 
to  hold  for  some  time  to  come. 
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FIGURE  CAPTIONS 


1.  Configuration  of  an  axial  vortex  in  an  axially  symmetric  strain  field. 

The  streamlines  of  the  strain  field  are  hyperbolic  in  each  cylindrical  plane. 
The  azimuthal  velocity  is  in  the  direction  of  the  curved  arrows. 

A 

2.  Swirl  velocity  profiles  of  various  times  t for  the  collapsing  strained 
vortex.  The  local  Reynolds  number  (e)”l  = 500. 

A 

3.  Mixing  configuration  at  various  times  t for  an  initially  planar  Interface 
displaced  one  half  vortex  radius  from  the  flow  symmetry  axis.  The  local 
Reynolds  number  (e)~^  = 500.  The  vortlclty  to  strain  rate  ratio  = S2/a  = 2. 

A 

4.  Mixing  configuration  at  various  times  t for  initially  circular  reactant 
displaced  one  half  radius  from  vortex  center.  Same  flow  parameters  as  Fig.  3. 

A A 

5.  Lagranglan  angular  coordinate  at  time  intervals  t of  0.5  from  t =0.5  to 

A 

t=5.0.  Same  flow  parameters  as  Fig.  3. 

6.  Sketch  of  mixing  and  diffusing  reactant  configuration  in  Lagranglan 
coordinates.  The  shaded  annulus  is  the  region  of  Intense  vortlclty,  the  Inner 
region  has  no  swirl  velocity,  and  the  outer  region  is  irrotatlonal . 
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SECTION  III 


A SOLUTION  FOR  A TWO-DIMENSIONAL  DIFFUSION  FLAME  IN  A VORTEX  FIELD 

The  theoretical  study  of  chemical  reactions  in  complex  flow  fields,  for 
example  turbulent  reacting  flows,  has  received  increased  attention  recently. 
This  attention  is  warranted  because  both  numerical  and  analytical  progress  is 
being  made  in  addressing  the  problem;  discussion  of  some  of  this  progress  is 
presented  in  Reference  1.  Also  presented  in  Reference  1 is  an  outline  of  a 
rather  general  approach,  using  both  analytical  and  numerical  methods,  for 
attacking  turbulent  reacting  flows  and  an  analysis  of  a model  of  small-scale 
mixing  and  reaction. 

In  Reference  2,  a similar  but  more  specialized  model  problem  of  small-scale 
mixing  and  reaction  has  been  posed  and  analyzed  locally  along  the  flame  front 
to  provide  dependencies  of  the  solution  upon  dimensionless  parameters.  This 
problem  has  several  important  features;  it  is  a diffusion-controlled  reaction 
in  a viscously  spreading  vorticity  field  which  stretches  the  flame  sheet.  The 
problem,  named  the  Marble  problem  here,  allows  rather  detailed  treatment,  as 
seen  in  Reference  2 which  has  slight  variations  in  detail  and  notation  from 
the  present  write-up. 

The  novel  features  of  this  work  compared  with  the  analysis  in  Reference  2 are 
(i ) it  is  observed  that  the  convection-diffusion  equation  for  the 
Shvab-Zel dovi ch  variable  permits  a similarity  solution,  reducing  the  number  of 
independent  variables  from  three  (radius,  angle  and  time)  to  two  (angle  and 
similarity  variable);  (ii)  as  in  Reference  1,  a Lagrangian  coordinate  system 
is  used  to  eliminate  flame-sheet  resolution  problems  induced  by  vortex 
winding;  (iii)  Fourier  analysis  in  angle  and  a numerical  treatment  in  the 
similarity  variable  allow  one  to  solve  the  global  problem.  The  analysis 
outlined  below  has  not  yet  been  completed.  The  solution  to  this  problem 
should  provide  a special  case  on  which  to  test  methods  for  solution  to  the 
more  general  problems  described  in  Referernce  1. 

Initially  we  have  fuel  in  the  left  half-plane  and  oxidizer  in  the  right 
half-plane  in  arbitrary  proportions  (see  Figure  1).  At  t=0,  the  two  half 
spaces  are  brought  into  contact  and  simultaneously  the  line  vortex  with 
tangential  velocity  ve  is  imposed: 


d0  r 

V0  (r,t)  = r --  = + -—  [l-exp(-n)]  (1) 

dt  2Trr 

where  r is  the  circulation  of  the  vortex,  v is  the  kinematic  viscosity  and 
n = r2/4vt.  (The  function  ve(r,t)  is  a solution  of  the  Navi er-Stokes 
equati ons ) . 

If  the  mass  fraction  of  fuel  in  the  left  half-space  initially  is  Yf  q while 
the  mass  fraction  of  oxidizer  in  the  right  half-space  initially  is  q,  and 
if  the  single-step  reaction  occurs 


vf[Yf]  + vo[Yq]  -»■  vp[Yp]  , 


then  with  Mf,  Mq  = the  molecular  weight  of  fuel  and  oxidizer  respectively, 


Z E { ----- 

vfMf 


^0  Yo,o  . Yf^o  Yq,o  . 

+ I / I + j 

\^o<^o  '^oMo  VfMf  vqMq 


is  the  Zhvab-Zel dovich  variable  for  which 


since 


Z = 1 for  X < 0 

Z = 0 for  X > 0 @ t = 0 


Yf  = Yf^o  > 

Yf  = 0 


Yq  = 0 for  X < 0 

Yo  = Yq,o  for  X > 0 @ t = 0. 


Z satisfies  the  convection-diffusion  equation: 


9Z 

--  + V*  (uZ)  = DV^Z 

9t 


or  in  cylindrical  coordinates 


9Z  V0(r,t)  9Z  92Z  1 9Z  1 92Z 

9t  r 90  9r2  r 9r  r2  302 

For  the  flame  sheet  approximation,  Yf  and  Yq  cannot  coexist,  so  that 


Yq,o  Yq,o  Yf^o 

Zf  = / { + } 

VqI^O  VqMq  VfMf 


determines  the  flame  surface. 

Integrating  the  tangenti al -vel ocity  equation  gives  the  angle  0(r,t,0o)at 
any  fluid  element  initially  located  at  r,0Q: 


(2) 


(3) 


(^) 


(5) 


(6) 


t for 


. r t 

- 00  = — [^  ~ ! (-r'2/4vt"  )dt" } 

ZttT^  0 


or 

r 1 

e(r,t)  - 00  = - [1  - E2(n)]  (7) 

2tt»4v  n 


where 


00  dt 

E2(z)  = / --  exp  (-zt) 

1 t2 


and 

Change  independent  variables, 
r = p 


n = r2/4vt 

to  Lagrangian  coordinates. 


r 1 

9 = 0Q  + - [1-E2(n)] 

2-rr»4v  n 

t = T 

Then  Eq.  (5)  becomes,  in  Lagrangian  coordinates. 


3Z  32z  1 3Z  2r  1 1 a2z 

..  = of—  + - - + f[—  - (1  - e-0)]2  t 1}  - — - 

9t  3p2  p 3p  8tt2  n p^  36q 

2r  2 8Z  4r  1 1 92z 

+ [e-^-(l-e-n)/n]  — + (1-e-n)  ---  } 

8irv  p2  39q  8ttv  p n 9p90o 


(8) 


(9) 


This  equation  admits  a similarity  solution  (which,  in  fact,  exactly  solves  the 
Marble  problem  since  neither  the  velocity  field  generated  by  the  line  vortex, 
nor  the  diffusion  field  produced  by  the  two  adjacent  half-spaces  of  fuel  and 
oxidizer  have  a length  scale  associated  with  them.) 


Let  fi  = p2/4vT  as  above,  and  assume  Z(p,t,9o)  = Z(n,0o). 
Then  Eq.  (9)  becomes 


32Z  3Z  2r  1 1 3^Z 

— - + (Sc  + 1/n)  — + {[ (1  - e-n)]^  + 1}  ---  — - 

3n2  3n  8irv  n 4n^  39q 


r 1 1 3Z  2r  1 3^Z 

+ [e-9  - - (1  - e-n)]  — - + — - — (1  - e-n)  = 0 (10) 

Sttv  n 39q  8itv  3h39q 


where  Sc  = v/D  is  the  Schmidt  number. 
Fourier  decomposition  of  Eq.  (10)  yields 


d^Zn  2r  1 dZn 
----  + [(Sc  + 1/n)  + in  — - --  (1  - e-n)]  --- 
dn^  Sttv  dn 


2r  1 


r 1 


+ {-{[  - (1  - e-n)]^  + 1}  ---  + i n [e-n  - - (1  - e-n)]}Zn  = 0 


8itv 


4n‘ 


8tt V n‘ 


(11) 


The  fact  that  a similarity  solution  exists  has  reduced  the  number  of 
independent  variables  upon  which  the  solution  depends  by  one.  When  Fourier 
analyzed,  the  mathematical  problem  is  reduced  from  a one-dimensional, 
time-dependent  partial  differential  equation,  an  i ni ti al -val ue  problem,  to  an 
ordinary  differential  equation  in  the  similarity  variable,  a two-point 
boundary-value  problem.  For  each  Zn(n)  we  have  a linear  ordinary  differential 
equation  for  which  boundary  conditions  are  needed  as  n-*-”  and  as  n>0.  For 
fixed  p,  as  t-^0,  n=P^/4vt-H»  and  Zn(n)  must  satisfy  initial  conditions  for  the  ^ 
Marble  problem: 


Z(n,eo) 


Zn(n)  exp(in0o) 


As  Z(r|x”,9o)  = H(9o“ir/2)  - H(9o-3tt/2) 

for  0 < 9g  < 2tt 


(12) 


= I Anexp(in9o) 

n=-oo 


(13' 


Hence 


Ao  = 1/2 
A2j  = 0, 


(-l)j 

A2i-1  = j = l,2,-3,4,... 

Tr(2j-1) 


A-2j  = 0 
A-(2j-l) 


(-l)j 

»(2j-l) 


(14) 


As  n>0,  Zn(n)  must  remain  finite.  Examination  of  the  indicial  equation  for 
Zn(n)  as  n>0  shows  that  the  two  solutions  behave  as  where 
u = -inr/STTv  + n/2.  Hence  a change  of  independent  variable  to 


Yn(n) 


nm/2 

l+nf^/2 


Zn(n) 


(15) 


where  m = ] n | , 

and  requiring  Yn(n)>0  as  n>0  implies  that  Zn(n)  remains  finite  at  n=0. 
The  equation  for  Yn(n)  is 


d^Yp  m nr  1 dYp 

— r-  + [(Sc  + 1/n) + i — - --  (1-e-n)]  — - 

dn“^  n(l+n'^/^)  4ttv  dn 


m/2(m/2+l)  m/2(Sc+l/n) 

n2(i+^m/2)  n(l+n'^/2) 


r 1 n^ 

j[  ...  . (i.e-'i)]2  + 1 i -- 

4ttv  n 


r 1 1 

+ i n — --  |[e"^  - - (l-e~'^)] 

8itv  n 


m 

l+n'"/2 


- (l-e-i)ll  Yn 


0 


(16) 


It  is  convenient  to  map  the  interval  0 < n < <»  into  0 < x < 1 using  the 
transformation 


X = 1-e"^ 


(17: 


Then  the  equation  for  Yn(x)  can  be  solved  using  boundary  conditions  that  Yq(x) 
= 0 at  X = 0 and  Yp(x)  = Ap  at  x = 1.  Equation  (16)  can  be  discretized 
strai ghtforwardly  using  a second-order,  central  difference  scheme,  and  the 
resulting  linear  algebraic  equations  can  be  solved  using  a tridiagonal  solver. 
The  Fourier  synthesis  to  obtain  the  complete  solution  can  be  performed  using 
an  FFT  package.  Location  of  the  flame  sheet  and  of  the  rate  of  fuel 
consumption  must  then  be  computed  from  the  complete  solution.  Work  is  in 
progress  on  these  tasks. 

Plots  of  the  interface,  0 = 0 and  9 = tt,  when  the  diffusion  0 = 0 for  four 
values  of  the  Reynold's  number  r/8irv  are  shown  in  Figures  2-5.  When  the 
diffusion  is  zero,  the  solution  depends  only  on  this  single  parameter.  It  is 
important  to  note  that  this  sequence  of  plots  does  not  represent  a time 
sequence;  rather  each  plot  shows  the  interface  for  that  Reynold's  number,  the 
similarity  variable  n as  a function  of  the  angle  0.  This  interface  is 
calculated  from  Eq.  (7). 
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